{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Set Algebra**\n",
    "\n",
    "$S \\cap \\left(T \\cup U\\right) = \\left(S \\cap T \\right) \\cup \\left(S \\cap U \\right)$\n",
    "\n",
    "$S \\cup \\left(T \\cap U \\right) = \\left(S \\cup T \\right) \\cap \\left(S \\cup U \\right)$\n",
    "\n",
    "**De Morgan's Laws**\n",
    "\n",
    "$$\\left( \\bigcup_n S_n \\right)^c = \\bigcap_n S_n^c$$\n",
    "\n",
    "Represents everything that is not in *at least one set* and in words reads as the compliment of the union of a collection of sets equals the intersection of the set compliments.\n",
    "\n",
    "$$\\left( \\bigcap_n S_n \\right)^c = \\bigcup_n S_n^c$$\n",
    "\n",
    "Represents everything that is not in all sets and in words reads as the compliment of the intersection of a collection of sets is equal to the union of the compliments."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Probability Axioms**\n",
    "\n",
    "1. **Nonnegativity** - For every event *A*, $P(A) \\ge 0$\n",
    "2. **Additivity** - For any sequence of **disjoint** events, the probability of their union is the sum of their individual probabilities <br/> <br/>\n",
    "$$P\\left( A_1 \\cup A_2 \\cup \\ldots \\right) = P\\left(A_1\\right) + P\\left(A_2\\right) + \\ldots$$ <br/>\n",
    "3. **Normalization** - The probability of the entire sample space is equal to 1"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Conditional Probability**\n",
    "\n",
    "$$P\\left(A \\vert B\\right) = \\frac{P\\left(A \\cap B\\right)}{P\\left(B\\right)}$$\n",
    "\n",
    "Can think of $P(B)\\hspace{1pt}$ as the fraction of the total sample space $\\Omega$ occupied by the event $B$ and $P(A\\cap B)\\hspace{1pt}$ as the fraction $\\Omega$ occupied by the interesection of $A$ and $B$, thus the conditinal probability is the fraction of $B\\hspace{1pt}$ occupied by the event $A\\cap B\\hspace{1pt}$. **NOTE** Conditional probabilities are themselves a valid probability law and hence all statements that hold for unconditional probabilities also hold for conditional probabilities."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Multiplication Rule**\n",
    "\n",
    "Assuming all conditioning events have positive probability, the following holds\n",
    "\n",
    "$$P\\left(\\cap_{i=1}^n A_i\\right) = P\\left(A_1\\right)P\\left(A_2 \\vert A_1 \\right) P\\left(A_3 \\vert A_1 \\cap A_2 \\right) \\ldots P\\left(A_n \\vert \\cap_{i=1}^{n-1} A_i \\right)$$\n",
    "\n",
    "Note there is nothing special about the ordering of the events in this statement. Any event, $A_i\\hspace{1 pt}$ can be used as the leading term see page 24 of text."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "*Example: Drawing three cards without a heart*\n",
    "\n",
    "This is example 1.10, page 25 from [1]. Consider three cards drawn from an ordinary 52-card deck without replacement. Find the probability that none of the three cards is a heart. As shown in the text this problem can be solved using the multiplication rule. Namely, define the event\n",
    "\n",
    "$$A_i = \\left\\{ \\mbox{the} \\hspace{1pt} i^{th} \\hspace{1pt} \\mbox{card is not a heart} \\right\\}$$\n",
    "\n",
    "The probability we seek is \n",
    "\n",
    "$$P\\left(A_1 \\cap A_2 \\cap A_3\\right) = P\\left(A_1\\right)P\\left(A_2 \\vert A_1 \\right) P\\left(A_3 \\vert A_1 \\cap A_2 \\right)$$\n",
    "\n",
    "As shown in the text this equal to \n",
    "\n",
    "$$P\\left(A_1\\right)P\\left(A_2 \\vert A_1 \\right) P\\left(A_3 \\vert A_1 \\cap A_2 \\right) = \\frac{39}{52}\\frac{38}{51}\\frac{37}{50}\\approx 0.414$$\n",
    "\n",
    "Below we simulate drawing three cards from a deck. The deck is represented as an array of integers from 1 to 52. Let integers 1-13 represent the cards with hearts. We \"shuffle\" the deck with the [Fisher and Yate's Algorithm](http://en.wikipedia.org/wiki/Fisher%E2%80%93Yates_shuffle) and examine the cummulative probability of the event of no hearts being drawn in the first three draws as well as the three events used in the calculation above."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[748, 540, 421]\n",
      "Theoretical Results:\n",
      "P(A1) = 0.75 \t P(A2 | A1) = 0.745098039216 \t P(A3 | A1,A2) = 0.74\n",
      "Observed Results:\n",
      "P(A1) = 0.748 \t P(A2 | A1) = 0.72192513369 \t P(A3 | A1,A2) = 0.77962962963\n",
      "Theoretical probability of no hearts in first three draws\n",
      "P(A1,A2,A3) = 0.413529411765\n",
      "Observed probability of no hearts in first three draws\n",
      "P(A1,A2,A3) = 0.421\n"
     ]
    },
    {
     "data": {
      "text/plain": [
       "[<matplotlib.lines.Line2D at 0x1048bd0d0>]"
      ]
     },
     "execution_count": 1,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": [
       "iVBORw0KGgoAAAANSUhEUgAAAXsAAAD9CAYAAABdoNd6AAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\n",
       "AAALEgAACxIB0t1+/AAAIABJREFUeJzt3XtcVGX+B/DPuEPlhbzibYYNBQQMRQyyixZdXLRy1DTD\n",
       "dm3XyPjtRlbbtrVtv033t5l02dVi26W2i9qqdMfcxGJtzPUCGmmtukoKNeAtSLymwPD9/fHdYRgH\n",
       "BtAZpjif9+s1L86Z88w5zzxz5nu+5znPHEwiIiAiog6tU7ArQEREgcdgT0RkAAz2REQGwGBPRGQA\n",
       "DPZERAbAYE9EZAA+g/0dd9yBfv36YdiwYc2WmT17NqKjo5GQkIBPP/3U7xUkIqJz5zPYz5w5E/n5\n",
       "+c0uf//99/HFF1+gpKQEL7zwAn7+85/7vYJERHTufAb7MWPGoGfPns0uX7FiBX76058CAEaNGoXq\n",
       "6mocPHjQvzUkIqJzZj6XF1dUVCA8PLxh3mq1ory8HP369fMoZzKZzmUzRESG5a+bHJzzBdozK9Jc\n",
       "YH/qKYEIH4899ljQ6/BdebAt2BZsC98PfzqnYG+xWOBwOBrmy8vLYbFYzrlSRETkX+cU7G02GxYv\n",
       "XgwA2LRpE3r06OHVhePi54MUERG1gc8+++nTp2Pt2rWorKxEeHg45s6di9raWgBARkYGbrjhBrz/\n",
       "/vuIiopC165d8corr7RLpb/PUlJSgl2F7wy2hRvbwo1tERgm8XfHUFMbMZmQlSX49a8DvSUioo7D\n",
       "ZDL5re+ev6AlIjKAdgv27LMnIgoeZvZERAbAzJ6IyACY2RMRGQAzeyIiA2BmT0RkAMzsiYgMgMGe\n",
       "iMgA2I1DRGQAzOyJiAyAmT0RkQEwsyciMgBm9kREBsDMnojIAJjZExEZADN7IiIDYGZPRGQAzOyJ\n",
       "iAyAmT0RkQEwsyciMgAGeyIiA2A3DhGRAbQY7PPz8xEbG4vo6GhkZWV5LT98+DAmT56MhIQEjBo1\n",
       "Ctu3b29yPczsiYiCx2ewdzqdyMzMRH5+Pnbs2IFly5Zh586dHmXmzZuHkSNHYtu2bVi8eDHuvffe\n",
       "gFaYiIjazuxrYVFREaKiohAREQEASEtLQ15eHuLi4hrK7Ny5Ew8//DAAICYmBmVlZfj6668RFhbm\n",
       "sa61a+dgzhydTklJQUpKit/eBBFRR2C322G32wOybp/BvqKiAuHh4Q3zVqsVhYWFHmUSEhLw9ttv\n",
       "Y/To0SgqKsKXX36J8vJyr2B/1VXuYE9ERN7OTITnzp3rt3X77MYxmUwtruDhhx9GdXU1EhMTkZ2d\n",
       "jcTERPzgBz/wKsc+eyKi4PGZ2VssFjgcjoZ5h8MBq9XqUSY0NBQvv/xyw/ygQYMwePBgP1eTiIjO\n",
       "hc/MPikpCSUlJSgrK0NNTQ1yc3Nhs9k8yhw5cgQ1NTUAgBdffBFXX301unXr5rUuZvZERMHjM7M3\n",
       "m83Izs5GamoqnE4n0tPTERcXh5ycHABARkYGduzYgZ/97GcwmUyIj4/HSy+91C4VJyKi1jOJBD7n\n",
       "NplMePhhwRNPBHpLREQdh8lkgr9CNH9BS0RkALw3DhGRATCzJyIyAGb2REQGwGBPRGQA7MYhIjIA\n",
       "ZvZERAbAzJ6IyACY2RMRGQAzeyIiA2BmT0RkAMzsiYgMgJk9EZEBMLMnIjIAZvZERAbAzJ6IyACY\n",
       "2RMRGQCDPRGRAbAbh4jIAJjZExEZADN7IiIDYGZPRGQALQb7/Px8xMbGIjo6GllZWV7LKysrMW7c\n",
       "OIwYMQLx8fF49dVXA1FPIiI6Bz6DvdPpRGZmJvLz87Fjxw4sW7YMO3fu9CiTnZ2NxMREbN26FXa7\n",
       "HQ888ADq6uq81sXMnogoeHwG+6KiIkRFRSEiIgIhISFIS0tDXl6eR5kBAwbg6NGjAICjR4+id+/e\n",
       "MJvNgasxERG1mc+oXFFRgfDw8IZ5q9WKwsJCjzKzZs3Ctddei4EDB+LYsWN4/fXXm1xXUdEczJmj\n",
       "0ykpKUhJSTmnihMRdTR2ux12uz0g6/YZ7E0mU4srmDdvHkaMGAG73Y49e/Zg7Nix2LZtG0JDQz3K\n",
       "JSe7gz0REXk7MxGeO3eu39btsxvHYrHA4XA0zDscDlitVo8yGzZswC233AIAiIyMxKBBg7Br1y6v\n",
       "dbHPnogoeHwG+6SkJJSUlKCsrAw1NTXIzc2FzWbzKBMbG4uCggIAwMGDB7Fr1y4MHjw4cDUmIqI2\n",
       "89mNYzabkZ2djdTUVDidTqSnpyMuLg45OTkAgIyMDDzyyCOYOXMmEhISUF9fjyeffBK9evXyWhcz\n",
       "eyKi4DGJBD4Mm0wmZGQI/vrXQG+JiKjjMJlM8FeI5u0SiIgMgLdLICIyAGb2REQGwMyeiMgAmNkT\n",
       "ERkAM3siIgNgZk9EZADM7ImIDICZPRGRATCzJyIyAGb2REQGwMyeiMgAGOyJiAyA3ThERAbAzJ6I\n",
       "yACY2RMRGQAzeyIiA2BmT0RkAMzsiYgMgJk9EZEBMLMnIjIAZvZERAbAzJ6IyABaDPb5+fmIjY1F\n",
       "dHQ0srKyvJY//fTTSExMRGJiIoYNGwaz2Yzq6uqAVJaIiM6OSaT5nNvpdCImJgYFBQWwWCxITk7G\n",
       "smXLEBcX12T5lStXYsGCBSgoKPDciMmE224T/P3v/q08EVFHZjKZ4CNEt4nPzL6oqAhRUVGIiIhA\n",
       "SEgI0tLSkJeX12z5pUuXYvr06U0uYzcOEVHwmH0trKioQHh4eMO81WpFYWFhk2VPnjyJ1atX4/nn\n",
       "n29y+eefz8GcOTqdkpKClJSUs6owEVFHZbfbYbfbA7Jun8HeZDK1ekXvvfceRo8ejR49ejS5PD7e\n",
       "HeyJiMjbmYnw3Llz/bZun904FosFDoejYd7hcMBqtTZZdvny5c124RARUXD5DPZJSUkoKSlBWVkZ\n",
       "ampqkJubC5vN5lXuyJEj+PjjjzFx4sRm18U+eyKi4PHZjWM2m5GdnY3U1FQ4nU6kp6cjLi4OOTk5\n",
       "AICMjAwAwLvvvovU1FR07tw58DUmIqI28zn00m8bMZkwbZogNzfQWyIi6jjabeglERF1DLxdAhGR\n",
       "ATCzJyIyAGb2REQGwMyeiMgAmNkTERkAgz0RkQGwG4eIyACY2RMRGQAzeyIiA2BmT0RkAMzsiYgM\n",
       "gJk9EZEBMLMnIjIAZvZERAbAzJ6IyACY2RMRGQAzeyIiA2BmT0RkAAz2REQGwG4cIiIDYGZPRGQA\n",
       "zOyJiAygxWCfn5+P2NhYREdHIysrq8kydrsdiYmJiI+PR0pKSpNlmNkTEQWP2ddCp9OJzMxMFBQU\n",
       "wGKxIDk5GTabDXFxcQ1lqqurcffdd2P16tWwWq2orKwMeKWJiKhtfGb2RUVFiIqKQkREBEJCQpCW\n",
       "loa8vDyPMkuXLsWUKVNgtVoBAH369GlyXczsiYiCx2dmX1FRgfDw8IZ5q9WKwsJCjzIlJSWora3F\n",
       "Nddcg2PHjuHee+/FjBkzvNb1xRdzMGeOTqekpDTb3UNEZFR2ux12uz0g6/YZ7E0mU4srqK2tRXFx\n",
       "Mf75z3/i5MmTuPzyy3HZZZchOjrao1xkpDvYExGRtzMT4blz5/pt3T6DvcVigcPhaJh3OBwN3TUu\n",
       "4eHh6NOnDzp37ozOnTvjqquuwrZt27yCPRERBY/PPvukpCSUlJSgrKwMNTU1yM3Nhc1m8ygzceJE\n",
       "/Otf/4LT6cTJkydRWFiIoUOHeq2LffZERMHjM7M3m83Izs5GamoqnE4n0tPTERcXh5ycHABARkYG\n",
       "YmNjMW7cOAwfPhydOnXCrFmzmgz2REQUPCaRwOfcJpMJY8cKPvgg0FsiIuo4TCYT/BWiebsEIiID\n",
       "4O0SiIgMgJk9EZEBMLMnIjIAZvZERAbAzJ6IyACY2RMRGQAzeyIiA2BmT0RkAMzsiYgMgJk9EZEB\n",
       "MLMnIjIAZvZERAbAzJ6IyACY2RMRGQCDPRGRAbAbh4jIAJjZExEZADN7IiIDYGZPRGQAzOyJiAyA\n",
       "mT0RkQEwsyciMoAWg31+fj5iY2MRHR2NrKwsr+V2ux3du3dHYmIiEhMT8Yc//KHJ9TCzJyIKHrOv\n",
       "hU6nE5mZmSgoKIDFYkFycjJsNhvi4uI8yl199dVYsWJFQCtKRERnz2dmX1RUhKioKERERCAkJARp\n",
       "aWnIy8vzKietSNuZ2RMRBY/PzL6iogLh4eEN81arFYWFhR5lTCYTNmzYgISEBFgsFjz99NMYOnSo\n",
       "17qKi+fghhuA5GTgmmtSkJKS4p93QETUQdjtdtjt9oCs22ewN5lMLa5g5MiRcDgc6NKlC1atWoVJ\n",
       "kyZh9+7dXuVqa+dg1Srg178GGOeJiLylpHgmwnPnzvXbun1241gsFjgcjoZ5h8MBq9XqUSY0NBRd\n",
       "unQBAIwfPx61tbX45ptvml1nff25VJeIiM6Gz2CflJSEkpISlJWVoaamBrm5ubDZbB5lDh482NBn\n",
       "X1RUBBFBr169ml0n++6JiNqfz24cs9mM7OxspKamwul0Ij09HXFxccjJyQEAZGRk4M0338Rf/vIX\n",
       "mM1mdOnSBcuXL/e5QWb2RETtzyStGUpzrhsxmQDoZvLzgdTUQG+RiOj7z2QytWq0Y2u0+y9o2Y1D\n",
       "RNT+gh7sT5wAli1r71oQERmLzz77QGjcZ19fD3TrptPTp7d3TYiIjKPdM/vGwf7yy9t760RExhS0\n",
       "bpx//hMoKmrvrRMRGVPQMvuvvvJ8vq6uvWtCRGQcQQn2x48Dp097Pr9nT3vXhIjIONr9Aq3TCVx8\n",
       "sXdmHxvLYZlERIHS7pl9XZ13oCcKli1bgCFDgAsvBBYvDnZtiAKn3YN9gO7eSdQmublATAxw003A\n",
       "yJFARgYwdy7w6KPa1cizTOpo2r0bZ+vW5pd9+y3QuXP71YW+W06d0rO+IUNaV/70aeD889u2jcOH\n",
       "gRtuALZtA554AoiM1IAP6O23J00CfvADnS8oAK67ruX1VVQA0dFtrwtRewr6L2h/+Uv39OHD7VsX\n",
       "8vbWW8C+fe756urA3rzuyy81wD79NNClCzBihAbZ5tTV6T60dSvQu7c+li4Fvv66+dccPw7MmQM8\n",
       "8IAG9m7ddNjvvfe6Az0AhIXpkOBly4BFi4CJE4HzztPXlZe7y+XnAyYTEB4O9OoFDBsGXHABMGuW\n",
       "bovoO0naAQDRr6jIJZdIwzQgsn69e3r79vaojfHU14uUlrqnhwzR9t6/37Pc55+LnHeeyBVXiOzd\n",
       "KzJ7tpa7/34Rp1PXUVMj8sYbInv2NL+tjAyRyZOb/jxPnhR58EGRp5/W9QwfrtswmUSmTRNZvVqk\n",
       "Tx+RDz/0fm15uciYMSKdOulrFi4UmTlTp3v2FHnySZFjx0R27xb54x9FqqpEfvpT9/41aJDIXXfp\n",
       "e2mNL78Uyc4WGTxYZMAAkY8/FvnHP0TCwkT+/GeR558X+fe/RU6dEtmwQesPiFx6qcizz4ocOSLi\n",
       "cHi2TWWlSFlZ67ZP5M8Q3e7BPjHR/eX78Y9dy/Xxr381/fr6epGtW9u+3W+/1S+cUR09KnL4sMjy\n",
       "5SJms8iCBSIvvqhtnZQkkpoqMnq0zs+eLXL11SLPPSdyww3uz+SWW0QuvNDzAA2IxMSIjBsn8te/\n",
       "urf3f/+nyyIj3eWuuEKDYX29lvnFLzRwupZPm6bL6urc61m7VgN+QYH7uY8/1vLjxon85jciS5a4\n",
       "l9XUiOTk6PLkZJELLtDp0FARq1Xk7bf1IHMuVq0S6dtXA/2mTU2Xqa/XA+avfiUSFaXb79JFE5yn\n",
       "nhKJiHC/75/8RNtm2jR9n0ePnlv9qGP6Xgf7ESPcO/zvfqfL33hD5997r+nXf/aZSP/+ns99/nnz\n",
       "XzoRDfKASLdu/nkPgfbnP4u8/vrZv/70aZHx4zXDrKrSNps61d3Wd93lnrbbNUBedpnOx8a6P5u6\n",
       "On39o4+6A/B//qPB6qOPRDIzRdatE7n3Xn19ly4aXP/2N11HQoKWF9GDTEKCHlTCwzVIDx4sUl2t\n",
       "2XVlpR4ImrJ2rQbWv/xFX9ezp8i773oeFM5UXy/y+9/rAa2+XuSdd0ROnDj7Nj1TWZnIzp2tK1tT\n",
       "o0H8s89EHnhA2/a3v9Wzqa+/1oPe44+7P6Mf/lDPTNavF1mzpvVnH9Sxfa+Dveu0HdCd3eW220QW\n",
       "L2769R9+qOVd2dnmze51NGfDhpbLnAunU2TOHD17OBf794tER2s9Bw7UTFxEg1p9vb73P/1Juy9O\n",
       "n/Z+/d693lm36zFggHZj/PKXWra8XOS119yvrawU+eornT56VOTgwbbV3ekUKSzUbgtAJD/fu8yp\n",
       "UyI2m8jdd4t0767lW2vtWvd7WbiwbXX7PqmvF8nL0wOh6/2ed552h1VWegf+TZtEsrJENm7UbqSO\n",
       "xHUGSMqfwT6od7284AL3dM+ezV+gPXRI/5aVAXFxQHJyy9v5z3/c0xUVgMXinv/iC73g1vgfn3/0\n",
       "kQ6/27nTPRpj506guFjvyPmnPwH33w90+u8l7eJiveg3Z46OIAkP99z+smU60uPSS/Xre/Ik0LWr\n",
       "dz2fegooKQESE4EBA7Qdzj9fR5pMmgRs367Le/YE7r4bCA3VNvz3v3X5Pffoen78Y+Cxx4A33gDM\n",
       "Zr3Y2KcPMGGCe1sWi5ZzcV3gBHS9oaEtt2tjnTrp+yss1PaOjfUuc/75QF6eTi9YoHVrrauu0guy\n",
       "rs+jozKZAJtNH65hn++8A/zjH/oZjhqln2PfvsDChbpP3Hgj8NBD+hlMmgT86Ee6X6elAUeO6EX2\n",
       "u+/W+QkT9KLynj26HwHAnXcC117bdH1EtE7+IgK8+y4wdKgOd62tBVau1PqcPq3vtU8foKpK92tX\n",
       "OxQW6nd93z7dn2fMAHbt0nWazXqRPixM31PPnvqaadPc38X6er14v22bPu9w6L504436+hMngL17\n",
       "gZAQoH9/oHt3/77v75p2D/a1te7ptgb70lIN9o2tXasjRiZO9Hx+507g8ceBTZuA9ev1w3ZZvBj4\n",
       "8ENg40b3c64df/16DTIA8OqrwJNPAj/5ic4nJADXX6/T776rX75Dh4C//U3HaAO6Y69fD9xxB3DZ\n",
       "ZXoQyc0Fbr8dmDcP+NWvtMyuXToSZfly4PPPgfh4fX82m+7wrm0kJmp9hw3TaUB3SBEdORIZCWze\n",
       "7P4SP/JI020YaE0F+jO1JdC7dPRAfyZXMjF1qj5+/3vg7bd1H/vySx29tnat+yB99KjuV4sW6ZDV\n",
       "ceM0eIWGAitWAK+/rkHyuus0cB49qvvOHXcA/frpQeLYMf1Ve1wc8NxzwHvvAePHayIyeLBn/ZxO\n",
       "YMMGTXaKioCBA4GLLtJEoqpK98f9+4FnnwWsVl13dbUG1oMHtR5vvqnlr7hCh1vffbd+93v21O/M\n",
       "W2/pQSA9Hdi9WwP6vn2aLISFaQJx6pQmGlVVetA4ckT3lccf1xFSPXtqPcLC9GA5ZYqOwgoJAe67\n",
       "T+PBCy/o9kWAmhodeXXZZVreYgEuuQS48kqNWeed524D1328zmZ/DqZ2/7eEF12kOy2gO3B6uk7/\n",
       "8Y+aIS9Y4P36Rx7RMdHZ2bqTdumiH0Z9vX4Ye/ZoQIyPd7/GZgN+9jN3Fv/ss/r8gQN6a+WvvtLh\n",
       "er166Q554YWaEdx8s7sOo0a578zZo4futD/6kQbdrCxgzRrNSIYPB0aP1n+3WFYGvPSSvsZq9Ryy\n",
       "B+iOPnWqe37UKD0gNVZbqzvlm29q/VwHot27db6wUA86H3ygf/nbBHL59lvdd1oKRCdP6ndm/XoN\n",
       "glu2aLJx6636HVu9GnjmGT0DGDsWWLdOg2lWlgbH4cP1u3fokH6fy8v1AHTggJ7BPvqoPh8Wpvvn\n",
       "rbdqwH/xReCaa1p3dn42qqo0+amq0oPZtde6D6AuGzdqxn/HHe4E6sQJPaBs2KDfeYdDz6D27tXl\n",
       "F1+s3//Bg3WZ671ZLJqI3XOP+5fY/uTPf0vY7sHelQ0DwGuvubsVXnlFf127aJH36++8U5dNmgT8\n",
       "z/9owN27V/8WFmq2AniO4Y+O1szmm2+A2bOBTz7R7KFXL10+fDjw29/qEX7NGuB//1eP9OPH6wd5\n",
       "/LhmLdu364Fi1CjPswNAs4GQEO9Tv8sv1/f2wgv65QD0y/D885qFAcDMmdoWV17p2dVC9F2xb59+\n",
       "R1avBpKSNNjPnq1nvs2dcdXX6/ehI3SHOJ0aP+rq3GfPBw7o9zY6Wr/T+/ZpL0J2tk5PmqTJ2kUX\n",
       "abK3cqX+diQxUZOzwYN1+RVX6IGlqgq46y6NGZ98omdiN9+svQgnTgBhYd/jYB8aqpk0oP3Lriz3\n",
       "nXc04K9Y4f16m01ft3cv8POfa3awbp2eji1c6C53/LgG4M6dNRM/dkx3vt699ZRuyxZ3lvzss9rn\n",
       "98wz7v7w557T7oi//13Xv2KFdsO4nDgB3HabnuaWl+vBANDtLF2qp5IffaRdM927685SVaU7B6AH\n",
       "m+Rk/dFQRITfmpeIgkxE48PSpXqg27tXDwzjxgE7dmhSeuON2gvRpYvGgORkDf45ORorTp7U+JKX\n",
       "p/HK6QROn/ZfsG/30Thms3vEQeOhlh99pCNOmjJqlI7GcL3u1lv1edfY6sZDOW+5xXsUzpVX6rC9\n",
       "Cy8USUnRETS7d+vol/x8LZuTo2Ufflhkxgx9Lj296focOnT2bcHRBkTUmNOpI6xcv7VwOkUqKnT0\n",
       "nT9DdFDuennnnTrd+AJtr156ytSUQ4e0G8XFatW/cXG6jlOn9GZWb76pZwtnuvJKPb06ehQYM0Zf\n",
       "ExWlrxs3Ts8WZs3SsjffDCxZotMPPth0fcLCWv9+z9QRTm+JyH86ddL45hoN16mT9ho0vijsl+20\n",
       "VCA/Px+xsbGIjo5GlqsDugmbN2+G2WzG22+/3eJGQ0O1z6/xjaP69AEqK93zs2a579Fy6JDnKBzX\n",
       "MMqRI4Hf/U7XM2KEni65LFvmnh4zxj09c6b+NZncQXvYMHcQTkrSC7XPPKPDxIiIOgKf1+ydTicy\n",
       "MzNRUFAAi8WC5ORk2Gw2xJ0x/tHpdOKhhx7CuHHjWtW/FBKi/eqNM3vXONv6ej2y5ecDkye7+74b\n",
       "jwF3BfuuXYHf/EanXVfVAeAPf9DxxS6jR2smX1jovkALaP/aZ5/psEgXk0mvCQwf3uLbICL63vAZ\n",
       "7IuKihAVFYWI/15NTEtLQ15enlewf+655zB16lRs3ry5dRs1a7BvnNmfd55euKiu1oB8+LBe0Hjl\n",
       "Fb3A6RpbHhPj7sZpLDFR13HypPdIgR499IdJZxo5Uh9nSk1t1dsgIvre8BnsKyoqEN7op6FWqxWF\n",
       "hYVeZfLy8rBmzRps3rz5vyNvmjKnYaq8PAWdO6d4ZPaAdqtUVmoWf+KEPjdzpg5zclm5Un+4caak\n",
       "JB0fb7Qf4RBRx2G322EP0H948hnsmw/cbvfddx/mz5/fMPi/+W6cORg5Un95Fx2tF2mjojxLhIW5\n",
       "f+jUWOOsvHHgb+z88/VCKxHR91VKSgpSGt3HZa7rp/l+4DPYWywWOByOhnmHwwHrGX0on3zyCdL+\n",
       "20FeWVmJVatWISQkBDabzWt9V12lwd5s9rxo6uIK9q5x6URE5B8+g31SUhJKSkpQVlaGgQMHIjc3\n",
       "F8saD3MBsNf1e2IAM2fOxIQJE5oM9IBemG3890yuYF9d7fk8L5YSEZ0bn8HebDYjOzsbqampcDqd\n",
       "SE9PR1xcHHJycgAAGRkZbdqYK8g3d9+OxsH+yiv1TnRvveV5wzIiImq7Fu/bNn78eIwfP97jueaC\n",
       "/CuvvOJ7Y2bPv2cKC9PbER8+rLf7XbBAL9R26dJSLYmIyJd2/QWtK8g3143Tp4/eN/7WW3W4pMUC\n",
       "rFrVfvUjIuqoghLsfWX2Lj16BL4+RERG0a7BvjUXaF1c/4yDiIjOXVCCfXOZfb9+7mlm9kRE/tOu\n",
       "wd51e4TmMvvG4+uZ2RMR+U9Qgn1zmX3j2ycwsyci8p92Dfau+zO35h/1MtgTEflPUIJ9c904gN7Z\n",
       "MjOz6ZudERHR2WlFju0/LXXjuDz3XODrQkRkJN+5zJ6IiPzvO9tnT0RE/hOUX9B2avd/c05EZGzt\n",
       "GnZd/wulvr49t0pEREHJsRnsiYjaF4M9EZEBMNgTERlAuwX7O+4AEhL0H40nJrbXVomICABMIiIB\n",
       "34jJhHbYDBFRh+LP2MlBkEREBsBgT0RkAAz2REQGwGBPRGQADPbtzG63B7sK3xlsCze2hRvbIjBa\n",
       "DPb5+fmIjY1FdHQ0srKyvJbn5eUhISEBiYmJuOSSS7BmzZqAVLSj4I7sxrZwY1u4sS0Cw+f9J51O\n",
       "JzIzM1FQUACLxYLk5GTYbDbExcU1lLn++usxceJEAMDnn3+OyZMn44svvghsrYmIqE18ZvZFRUWI\n",
       "iopCREQEQkJCkJaWhry8PI8yXbt2bZg+fvw4+vTpE5iaEhHR2RMf3njjDbnzzjsb5pcsWSKZmZle\n",
       "5d555x2JjY2V7t27S2FhoddyAHzwwQcffJzFw198duOYXPckbsGkSZMwadIkrFu3DjNmzMCuXbs8\n",
       "lgt/PUtEFFQ+u3EsFgscDkfDvMPhgNVqbbb8mDFjUFdXh6qqKv/VkIiIzpnPYJ+UlISSkhKUlZWh\n",
       "pqYGubm5sNlsHmX27NnTkLkXFxcDAHr37h2g6hIR0dnw2Y1jNpuRnZ2N1NRUOJ1OpKenIy4uDjk5\n",
       "OQCAjIwMvPXWW1i8eDFCQkLQrVs3LF++vF0qTkREbeC33v9mrFq1SmJiYiQqKkrmz58f6M0F3Vdf\n",
       "fSUpKSkydOhQufjii2XhwoUiIlJVVSXXX3+9REdHy9ixY+Xw4cMNr5k3b55ERUVJTEyMrF69OlhV\n",
       "D4i6ujoZMWKE3HTTTSJi3HYQETl8+LBMmTJFYmNjJS4uTjZt2mTI9pg3b54MHTpU4uPjZfr06XLq\n",
       "1CnDtMPMmTOlb9++Eh8f3/Dc2bz3LVu2SHx8vERFRcns2bNbte2ABvu6ujqJjIyU0tJSqampkYSE\n",
       "BNmxY0cgNxl0+/fvl08//VRERI4dOyZDhgyRHTt2yIMPPihZWVkiIjJ//nx56KGHRERk+/btkpCQ\n",
       "IDU1NVJaWiqRkZHidDqDVn9/e+aZZ+S2226TCRMmiIgYth1ERG6//XZ56aWXRESktrZWqqurDdce\n",
       "paWlMmiqqvX5AAADxUlEQVTQIDl16pSIiEybNk1effVVw7TDxx9/LMXFxR7Bvi3vvb6+XkREkpOT\n",
       "G0Y+jh8/XlatWtXitgMa7Dds2CCpqakN80888YQ88cQTgdzkd87EiRPlww8/lJiYGDlw4ICI6AEh\n",
       "JiZGRPTI3fiMJzU1VTZu3BiUuvqbw+GQ6667TtasWdOQ2RuxHUREqqurZdCgQV7PG609qqqqZMiQ\n",
       "IfLNN99IbW2t3HTTTfLBBx8Yqh1KS0s9gn1b3/u+ffskNja24flly5ZJRkZGi9sN6L1xKioqEB4e\n",
       "3jBvtVpRUVERyE1+p5SVleHTTz/FqFGjcPDgQfTr1w8A0K9fPxw8eBAAsG/fPo8RTh2pje6//348\n",
       "9dRT6NTJvZsZsR0AoLS0FGFhYZg5cyZGjhyJWbNm4cSJE4Zrj169euGBBx7AD3/4QwwcOBA9evTA\n",
       "2LFjDdcOjbX1vZ/5vMViaVWbBDTYt3acfkd0/PhxTJkyBQsXLkRoaKjHMpPJ5LNtOkK7rVy5En37\n",
       "9kViYmKzv7MwQju41NXVobi4GL/4xS9QXFyMrl27Yv78+R5ljNAee/bswYIFC1BWVoZ9+/bh+PHj\n",
       "eO211zzKGKEdmtPSez8XAQ32bR2n31HU1tZiypQpmDFjBiZNmgRAj9gHDhwAAOzfvx99+/YF4N1G\n",
       "5eXlsFgs7V9pP9uwYQNWrFiBQYMGYfr06VizZg1mzJhhuHZwsVqtsFqtSE5OBgBMnToVxcXF6N+/\n",
       "v6HaY8uWLbjiiivQu3dvmM1m3Hzzzdi4caPh2qGxtnwnrFYrLBYLysvLPZ5vTZsENNi3Zpx+RyMi\n",
       "SE9Px9ChQ3Hfffc1PG+z2bBo0SIAwKJFixoOAjabDcuXL0dNTQ1KS0tRUlKCSy+9NCh196d58+bB\n",
       "4XCgtLQUy5cvx7XXXoslS5YYrh1c+vfvj/DwcOzevRsAUFBQgIsvvhgTJkwwVHvExsZi06ZN+Pbb\n",
       "byEiKCgowNChQw3XDo219TvRv39/XHjhhSgsLISIYMmSJQ2v8ckfFxx8ef/992XIkCESGRkp8+bN\n",
       "C/Tmgm7dunViMpkkISFBRowYISNGjJBVq1ZJVVWVXHfddU0Or3r88cclMjJSYmJiJD8/P4i1Dwy7\n",
       "3d4wGsfI7bB161ZJSkqS4cOHy+TJk6W6utqQ7ZGVldUw9PL222+Xmpoaw7RDWlqaDBgwQEJCQsRq\n",
       "tcrLL798Vu/dNfQyMjJS7rnnnlZt2yTCG9cQEXV0/E9VREQGwGBPRGQADPZERAbAYE9EZAAM9kRE\n",
       "BsBgT0RkAP8PL67/vVxmNsoAAAAASUVORK5CYII=\n"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "from functools import partial\n",
    "def swap(i,deck):\n",
    "    j = numpy.random.randint(0,i)\n",
    "    tmp = deck[j]\n",
    "    deck[j] = deck[i]\n",
    "    deck[i] = tmp\n",
    "    return deck\n",
    "\n",
    "def shuffledDeck(size=52):\n",
    "    deck = range(1,size+1)\n",
    "    idx = [size - i for i in range(1,size)]\n",
    "    map(partial(swap, deck=deck), idx)\n",
    "    return deck\n",
    "    \n",
    "totalExperiments = 1000\n",
    "\n",
    "def experiment(A):\n",
    "    deck = shuffledDeck()\n",
    "    hearts = range(1,14)\n",
    "    a1 = deck[0] not in hearts\n",
    "    a2 = deck[1] not in hearts and a1\n",
    "    a3 = deck[2] not in hearts and a1 and a2\n",
    "    return [a1+A[0],a2+A[1],a3+A[2]]\n",
    "\n",
    "rslt = [0,0,0]\n",
    "cumA = []\n",
    "for i in range(totalExperiments):\n",
    "    rslt = experiment(rslt)\n",
    "    cumA.append(rslt[2]/float(i+1))\n",
    "\n",
    "print rslt \n",
    "print \"Theoretical Results:\"\n",
    "print \"P(A1) = {0} \\t P(A2 | A1) = {1} \\t P(A3 | A1,A2) = {2}\".format(39/float(52),38/float(51),37/float(50))\n",
    "print \"Observed Results:\"\n",
    "print \"P(A1) = {0} \\t P(A2 | A1) = {1} \\t P(A3 | A1,A2) = {2}\".format(rslt[0]/float(totalExperiments), rslt[1]/float(rslt[0]), rslt[2]/float(rslt[1]))\n",
    "print \"Theoretical probability of no hearts in first three draws\"\n",
    "print \"P(A1,A2,A3) = {0}\".format(39/float(52)*38/float(51)*37/float(50))\n",
    "print \"Observed probability of no hearts in first three draws\"\n",
    "print \"P(A1,A2,A3) = {0}\".format(rslt[2]/float(totalExperiments))\n",
    "pyplot.plot(cumA)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Total Probability Theorem**\n",
    "Let $A_1, \\ldots, A_n\\hspace{1 pt}$ be disjoint events that form a partition of the sample space, i.e. each possible event is included in exactly one of the events $A_1, \\ldots, A_n\\hspace{1 pt}$, and assume $P(A_i)>0\\hspace{1 pt}$ for all *i*, then for any event *B*\n",
    "\n",
    "$$\\begin{eqnarray} P\\left( B \\right) &=& P\\left(A_1\\cap B\\right) + \\ldots + P\\left(A_n \\cap B\\right) \\cr\n",
    "                                &=& P\\left(A_1\\right)P\\left(B \\vert A_1\\right) + \\ldots + P\\left(A_n\\right)P\\left(B\\vert A_n\\right)\n",
    "\\end{eqnarray}$$\n",
    "\n",
    "A simple consequence of this theorem is that for a given set *A* we have the following\n",
    "\n",
    "$$P\\left(B\\right) = P\\left(A \\cap B \\right) + P\\left(A^c \\cap B \\right)$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Bayes' Rule**\n",
    "\n",
    "Let $A_1, A_2, \\ldots, A_n\\hspace{1 pt}$ be **disjoint** events that form a partition of the sample space and assume that $P(A_i)>0\\hspace{1 pt}$ for all *i*, then for any event *B* such that $P(B)>0\\hspace{1 pt}$, we have\n",
    "\n",
    "$$\\begin{eqnarray}\n",
    "P\\left(A_i \\vert B \\right) &=& \\frac{P\\left(A_i\\right) P\\left(B \\vert A_i\\right)}{P\\left(B\\right)} \\cr\n",
    "&=& \\frac{P\\left(A_i\\right) P\\left(B \\vert A_i\\right)}{P\\left(A_1\\right)P\\left(B\\vert A_1\\right) + \\ldots + P\\left(A_n\\right)P\\left(B\\vert A_n\\right)}\n",
    "\\end{eqnarray}$$\n",
    "\n",
    "where the second equality is a consequence of the Total Probability Theorem. \n",
    "\n",
    "Bayes' rule is often used for **inference**. There are a number of \"causes\" that may result in a certain \"effect\". We observe the effect and wish to infer the cause. Given that the effect *B* has occurred, we wish to evaluate the probability $P(A_i\\vert B)\\hspace{1 pt}$ that the cause $A_i\\hspace{1 pt}$ is present. Refer to $P(A_i \\vert B)\\hspace{1 pt}$ as the **posterior probability** of the event $A_i\\hspace{1 pt}$ and $P(A_i)\\hspace{1 pt}$ as the **prior probability**."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "*Example: The False-Positive Puzzle*\n",
    "\n",
    "This is example 1.18, page 33 from [1].\n",
    "\n",
    "Assume a test for a certain *rare* disease is 95% accurte in both postive and negative cases. Assume also that it is known that the 0.1% of the general population has the disease. This requires some clarification. The proper interpretation here is that, given the person has the disease, the test result is positve 95% of the time and given that the person does not have the disease, the test is negative 95% of the time. As we will see, this is **not** the same as saying that the test is 95% accurate accross all patients. Given that a person tests postive for the disease, what is the probability of having the disease?\n",
    "\n",
    "As shown in the text this can be solved using Baye's Rule as follows: Let *A* be the event that the person has the disease, and *B* be the event the test result is positive. The solution is then\n",
    "\n",
    "$$P\\left(A \\vert B \\right) = \\frac{P\\left(A\\right)P\\left(B\\vert A\\right)}{P\\left(A\\right)P\\left(B\\vert A\\right) + P\\left(A^c\\right)P\\left(B \\vert A^c\\right)}\n",
    "=\\frac{0.001 \\cdot 0.95}{0.001\\cdot 0.95 + 0.99 \\cdot 0.05}=0.0187$$\n",
    "\n",
    "Below we simulate this result:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Observed Values\n",
      "P(A) = 0.001 \t P(B) = 0.0493 \t P(A|B) = 0.0141987829615\n"
     ]
    }
   ],
   "source": [
    "numberPatients = 10000\n",
    "PA = 0.001 #probability that a random patient has the disease\n",
    "PBP = .95 #probability that test is positive if person actually has the disease\n",
    "\n",
    "ACount = 0 #number of times person actually has disease\n",
    "BCount = 0 #number of times test is positive\n",
    "ABCount = 0 #number of times person has disease A test is positive\n",
    "\n",
    "for i in range(numberPatients):\n",
    "    r = numpy.random.uniform(0,1,2)\n",
    "    if r[0]<=PA: #person has disease\n",
    "        ACount += 1\n",
    "        if r[1]<=PBP: \n",
    "            ABCount += 1\n",
    "            BCount += 1\n",
    "    else: #person doesn't have disease\n",
    "        if r[1]>PBP: BCount += 1 \n",
    "            \n",
    "print \"Observed Values\"\n",
    "print \"P(A) = {0} \\t P(B) = {1} \\t P(A|B) = {2}\".format(ACount/float(numberPatients),BCount/float(numberPatients),ABCount/float(BCount))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Independence**\n",
    "\n",
    "Two events *A* and *B* are defined as **independent** iff\n",
    "\n",
    "$$P\\left(A\\cap B \\right) = P\\left(A\\right)P\\left(B\\right)$$\n",
    "\n",
    "The above definition is equivalent to defining the events as independent when \n",
    "\n",
    "$$P\\left(A \\vert B\\right) = P\\left(A\\right)$$\n",
    "\n",
    "We use the first definition as it holds even when $P(B)=0 \\hspace{1 pt}$. Note that *set disjointness* does **not** imply independence. Indeed the opposite is true. Given disjoint sets *A* and *B* each with probability greater than 0, then the occurrance or non-oncurrance of either event provides complete information about the occurrance of the other event, e.g. $P(A|B) = 0\\hspace{1 pt}$, so that\n",
    "\n",
    "$$P\\left(A\\cap B \\right) = 0 \\lt P(A)P(B)$$\n",
    "\n",
    "**Several** events $A_1, A_2, \\ldots, A_n \\hspace{1pt}$ are independent if\n",
    "\n",
    "$$P\\left(\\bigcap_{i\\in S} A_i \\right) = \\prod_{i\\in S} P\\left(A_i\\right)$$\n",
    "\n",
    "for every subset $s \\subset S = \\left\\{1,2,\\ldots,n\\right\\} \\hspace{1 pt}$ Note it is **not** enough that the condition hold for *S*, it must hold for every subset of S. \n",
    "\n",
    "**Conditional Independence**\n",
    "\n",
    "Given an event *C*, the events *A* and *B* are called **conditionally independent** if\n",
    "\n",
    "$$P\\left(A \\cap B \\vert C \\right) = P\\left(A \\vert C\\right) P\\left(B \\hspace{1pt}\\vert\\hspace{1pt} C\\right)$$\n",
    "\n",
    "The above definition is equivalent to \n",
    "\n",
    "$$P\\left(A \\cap B \\vert C\\right) = P\\left(B\\vert C\\right) P\\left(A\\vert B \\cap C\\right)$$\n",
    "\n",
    "which can be interpreted as indicating that if *C* is known to have occurred, the additional information that *B* has occurred does not change the probability of *A*. Note, unconditional independence of *A* and *B* does **not** imply conditional independence and vice versa."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 2",
   "language": "python",
   "name": "python2"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 2
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython2",
   "version": "2.7.10"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 0
}
